Topological liquid nucleation induced by vortex-vortex interactions in Kitaev's honeycomb model 
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We provide a comprehensive microscopic understanding of the nucleation of topological quantum liquids, a 
general mechanism where interactions between non-Abelian anyons cause a transition to another topological 
phase, which we study in the context of Kitaev's honeycomb lattice model. For non-Abelian vortex excita- 
tions arranged on superlattices, we observe the nucleation of several distinct Abelian topological phases whose 
character is found to depend on microscopic parameters such as the superlattice spacing or the spin exchange 
couplings. By reformulating the interacting vortex superlattice in terms of an effective model of Majorana 
fermion zero modes, we show that the nature of the collective many-anyon state can be fully traced back to the 
microscopic pairwise vortex interactions. Due to RKKY-type sign oscillations in the interactions, we find that 
longer-range interactions beyond nearest neighbor can influence the collective state and thus need to be included 
for a comprehensive picture. The omnipresence of such interactions implies that corresponding results should 
hold for vortices forming an Abrikosov lattice in a p-wave superconductor, quasiholes forming a Wigner crystal 
in non-Abelian quantum Hall states or topological nanowires arranged in regular arrays. 

PACS numbers: 05.30.Pr, 75.10.Jm 



I. INTRODUCTION 

One of the most intriguing aspects of a topological phase 
is the emergence of anyonic quasiparticles. If these obey 
non-Abelian statistics, their presence gives rise to a (macro- 
scopic) ground state degeneracy, a distinctive feature which 
has been suggested to be exploited for topological quantum 
computation.' However, in any microscopic system this de- 
generacy will be lifted by the omnipresent interactions be- 
tween the anyons. These interactions are often assumed to 
be extremely weak due to the localized nature of the anyon 
wavefunction, but their strength grows exponentially when the 
anyons are brought into proximity. In fact, they can reach 
sizable magnitude when the anyon separation becomes of the 
order of a characteristic length scale, such as the magnetic 
length in quantum Hall liquids,^ the coherence length in a p- 
wave superconductor^'*, or the plaquette spacing in Kitaev's 
honeycomb model.!^ 

When interacting anyons form regular arrangements (e.g. 
a Wigner crystal in a fractional quantum Hall state'^ or an 
Abrikosov lattice in a topological superconductor), it has been 
shown that the collective degeneracy is split, and another topo- 
logical state (distinct from the parent state of which the anyons 
are excitations) is selected as the new ground state. ^ This 
mechanism through which local microscopies can change the 
global topological properties is referred to as topological liq- 
uid nucleation. Here we study it in the context of a micro- 
scopic model and show that the character of the nucleated 
topological state can be fully traced back to the signs and the 
relative magnitudes of the anyon-anyon interactions inherent 
in the model. This direct connection between pairwise inter- 
actions and the collective many-anyon state is made explicit 
through an effective Majorana zero mode lattice model for 
the interacting anyon lattice. Our results provide a compre- 



hensive understanding how nucleation occurs through a hy- 
bridization of localized modes and what states can emerge 
for different anyon lattice densities. Due to the general na- 
ture of the employed effective model, our results apply also 
to microscopically distinct systems supporting localized Ma- 
jorana mode such as fractional quantum Hal l liq uids^j p-wave 
superconductors'^ and topological nanowired^l^. 



II. VORTEX LATTICES IN KITAEV'S HONEYCOMB 
MODEL 

As a prototypical system that allows for good control of all 
microscopic parameters, we study nucleation in the context of 
Kitaev's honeycomb lattice model.'^It is an exactly solvable 
spin model defined by the Hamiltonian 

H= J, a] a] + kY. ' (1) 

7 — links {i,j,k) 

where the denote the standard Pauli matrices describing 
spin-1/2 moments on the sites of the lattice, 7 = x,y, z in- 
dicates the three different link types, are the strengths of 
the nearest neighbor spin exchange along these links and K is 
the magnitude of a three spin term that explicitly breaks time 
reversal symmetry. For a system of 2N spins the Hamiltonian 
([T]i has N local symmetries [H, Wp] = 0, where Wp are mutu- 
ally commuting Z2 valued six-spin operators associated with 
every plaquette p. Their eigenvalues Wp = —1 denote a vr-flux 
vortex at plaquette p, while Wp — 1 denotes an absence of 
one. One can thus restrict to a particular vortex sector labeled 
by the pattern of the eigenvalues {wp}. In each sector the 
spin model ([T]i can be mapped to a tight-binding model of free 
Majorana fermions tunneling on the honeycomb lattice. The 



corresponding Hamiltonian H^.^^j will always be quadratic in 
the Majorana operators and thus readily diagonalized.'^ 

The ground state over all vortex sectors resides in the 
vortex-free sector (Wp = 1 on all plaquettes), which supports 
phases with both Abelian Toric Code anyons as well as non- 
Abelian Ising anyons." We are interested here in the collec- 
tive properties of the latter, which appear as magnetic vortex 
excitations of a phase occurring in the vicinity of isotropic 
spin exchange Jx — Jy ^ Jz and for a finite three-spin ex- 
change > 0. By studying sectors with only two vortices, 
it has been shown that the vortices can combine into two dis- 
tinct collective states whose energies, due to a microscopic 
interaction, sensitively depend on the separation d between 
the vortices.^ The key properties of this pairwise interaction 
are summarized in Fig. [T] which illustrates how the energy 
splitting e((i) between the two collective states decays expo- 
nentially with increasing separation, i.e. e{d) ^ e^'^/^ with 
^ being the coherence length of the non-Abelian phase. It 
further shows RKKY-like sign oscillations at the wavelength 
of the inverse Fermi momentum. These characteristic fea- 
tures are not unique to the honeycomb model, but they have 
also been found for Ising anyons emerging in the Moore-Read 
quantum Hall state,'^^^; + ipy superconductors,^ and topolog- 
ical nanowires.'^The precise values for the amplitude and the 
frequency of the oscillations, however, depend always on the 
specific microscopic situation. 

The oscillating energy splitting e{d) is a property of a two- 
vortex problem. Our aim is to use it as an input to under- 
stand the many-vortex problem and answer which collective 
ground state is formed as we arrange the interacting vortices 
in superlattices. In the honeycomb model this corresponds 
to studying vortex sectors where the eigenvalues Wp = —1 
form a periodic pattern. For simplicity we restrict to con- 
sidering uniform and isotropic superlattices, which, as illus- 
trated in Fig.|4] we can parametrize with the superlattice spac- 
ing D = 1, 2, 3, . . . (in units of the plaquette spacing). All 
these sectors are translationally invariant with respect to a suit- 
ably chosen magnetic unit cell, with the corresponding Bloch 
Hamiltonians i/ij^iu,^} being readily (numerically) diagonal- 
izable 4£)^ x AD"^ matrices.E^I 



III. NUCLEATED PHASES AND VORTEX BAND 
STRUCTURE 

Varying the spacing D is of interest for two physical rea- 
sons. First, it enables us to tune the microsccopic pairwise 
interactions, as captured by the energy splittings e{D), which 
we find resulting in different collective topological states. We 
characterize them by the Chern number' ' that can be numer- 
ically evaluated for the ground state of each sector.EI Sec- 
ond, going systematically through the sequence of possible 
superlattice spacings enables us to simulate the expansion of a 
Wigner crystal of quasiholes when shifting the magnetic field 
on a quantum Hall plateau, or similarly the expansion of an 
Abrikosov lattice of vortices in a p-wave superconductor when 
varying the applied magnetic field. 

Our findings are given in Fig. |2] which shows two gen- 
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FIG. 1. Microscopies of the interaction between a pair of vortices 
in Kitaev's honeycomb lattice model: In agreement with the Ising 
anyon fusion rule ct x cr = 1 + -(/), the non-Abelian vortices can 
combine into two possible collective states with energies Ei and 
E^. The splitting e{d) = E^ — Ei between them decreases ex- 
ponentially with vortex separation and shows characteristic oscilla- 
tion due to interference effects. When e{d) > (< 0) the micro- 
scopies of the system energetically favours the state where the vor- 
tices combine into the trivial (fermionic) excitation. The plot is for 
Jx = Jy ~ Jz = 1 and K — 1/20, which corresponds to a coher- 
ence length of ^ ~ 2.5. The continuous curve has been obtained by 
suitably tuning the spin exchange couplings to simulate the continu- 
ous transport of the vorticesP 



eral ways the presence of a vortex superlattice modifies the 
phase diagram. First, the non-Abelian phase (characterized 
by Chern number = — 1) occurring in the vicinity of the 
isotropic spin exchange = Jy = Jz is always replaced by 
a distinct Abelian topological phase (characterized by an even 
Chern number). In particular, we find phases characterized by 
1/ — —2, —4, for integer superlattice spacings D = 1,2, 3, 
respectively. Second, we find that the Abelian phases existing 
in the dimerized limits (e.g. Jz > Jx + Jy), that are char- 
acterized hy ly = and support non-chiral Abelian anyons 
(so-called semionsjSE3 are always enlarged in the presence 
of a vortex superlattice. We focus here to understand the first 
phenomenon and show that all the Abelian phases emerging 
in the vicinity of Jx = Jy = Jz are nucleated phases, i.e. 
they arise due to the vortex-vortex interactions. 

A further investigation of the nucleated phases reveals char- 
acteristic band structures. The spectrum of the v — —2 phase 
arising in the presence of a full-vortex {D = 1) superlattice, 
shown in Fig. pi exhibits four gapped Dirac cones.IElThe spec- 
trum of the V = —4 phase arising in the presence of a Z) = 2 
superlattice also exhibits four gapped Fermi points, but with a 
significantly broader (quadratic) dispersion. The energy spec- 
trum of the nucleated phase for the D — 'i superlattice, how- 
ever, does not show indications of Fermi points in agreement 
with the phase being characterized hy v = Q. Unlike the 
other nucleated phases, this phase remains gapped even as the 
three-spin coupling K is tuned to zero. This indicates that it 
is adiabatically connected to the previously observed gapped 
phase,^^ that emerges when a vortex superlattice is imposed 
on the gapless time-reversal symmetric spin liquid (for which 
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FIG. 2. Phase diagram of Kitaev's honeycomb model (a) in the ab- 
sence of a vortex superlattice and (b) in the presence of D = 1, 2, 3 
superlattices. We have parametrized the phase diagram by J]] J.^ = 
1 and set K = 0.03. The plots show the lowest energy gap and 
Chem numbers v along a cut (1/3 < Jz < 1/2 when J^, = Jy) 
indicated by the dashed red line in (a). In the absence of vortices the 
continuous phase transition between the Abelian toric code phase 
(i^ = 0, denoted by TC) and the non- Abelian Ising phase {u = —1) 
occurs at Jz = 1 /2. In the presence of a vortex superlattice we again 
find a continuous phase transition, but which occurs now always at 
some Jf < 1/2 separating the TC phase (Jz > Jf ) from the nucle- 
ated Abelian phase (Jz < Jz). The exact value of Jz{D) can also 
be traced back to the vortex-vortex interactions, with there being no 
simple direct relation between the spacing D and the critical value.^" 



K = 0). Despite this behavior suggesting a different ori- 
gin (the non- Abelian anyons underlying the nucleated phases 
emerge only for K ^ 0), we will show below that also this 
phase can be traced back to the vortex-vortex interactions. 

The crucial common feature of all these observed band 
structures is that they consist of a low-energy vortex band ^ J 
and a set of high-energy /erm/on bands (see Fig.jsjfor an 
illustration). In the presence of 2N vortices the first contains 
N states that have support only on the sites near the vortices, 
while the latter contains the rest of the states that in general 
have support on all sites of the honeycomb lattice. For any 
finite K the vortex and fermion bands are separated in energy 
by a band gap A^f. Hence the Chem number characterizing 
the ground state can be written as 



(2) 



where and Vf are the Chern numbers for the occupied neg- 
ative energy bands and bands, respectively. We ob- 
serve that the first depends on the underlying vortex configura- 
tion, while the latter contributes always — ~1 for K > 0. 
The nucleated phases can thus be viewed as comprising of two 
decoupled theoiies: a remnant of the non- Abelian phase liv- 
ing on the honeycomb lattice and an emergent theory living 
effectively "on top of it" on the vortex lattice. The problem of 
understanding the nature of the collective many-vortex state 
is therefore reduced to the problem of understanding how 
depends on D. 
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FIG. 3. The characteristic band structure of the nucleated topo- 
logical phases with (*^/ ) denoting the low-energy (high-energy) 
vortex (fermion) bands. The actual plot is for the D = 1 vortex su- 
perlattice on a cylinder (open boundary conditions in x-direction). 
Consistent with the nucleated Chem number u — —2 phase, the 
spectral flow shows two edge modes per edge (red and green lines de- 
noting different edges) emanating from the four gapped Dirac cones. 



IV. EFFECTIVE MAJORANA MODEL FOR THE 
INTERACTING VORTEX SUPERLATTICE 

We now turn to address the microscopic connection be- 
tween the presence of a given vortex superlattice and the na- 
ture of the nucleated phase. To this end we consider an ef- 
fective tight-binding model of Majorana fermions that con- 
nects the energy splitting e{d) due to vortex-vortex interac- 
tions to the Chern number of the vortex band character- 
izing the nucleated topological phase. The motivation for 
such a model is as follows: The vortices bind localized Majo- 
rana zero modes and the oscillating energy splitting shown in 
Fig.[T]can be viewed as arising from them tunneling between 
the two vortex cores. ^ Likewise, we assume that the voitex 
band can be understood as arising through a collective 
hybridization, which occurs as the Majorana zero modes start 
tunnelling on the triangular vortex superlattice. To capture 
also the case when longer range interactions become relevant 
our model incorporates both nearest (n.n.) and next-nearest 
(n.n.n.) neighbor hopping of amplitudes ti and t^, respec- 
tively (see Fig.|4|. The effective Majorana Hamiltonian is then 
given by 



Hm ~ it 



n.n. 



2^ S^3 1^1■J 



it 



n.n.n. 



(3) 



where the 7^ denote Majorana zero modes at vortex location i 
(the center of honeycomb plaquette) obeying {7i , 7j } = 2Sij . 
The model allows for Z2 gauge choices Sij = ±1, which 
give rise to flux = —i\n{isijSjkSki) = ±7r/2 through 
each triangular plaquette with corners k. As illustrated 
in Fig. |4] there are three distinct types of plaquettes: those 
spanned solely by either ti- or i^g-links and those consisting 
of both. We denote the former as Ti and T^, respectively, 
and the latter as When = 0, it has been shown that 

for a uniform triangular lattice of 7r-flux vortices one should 
impose a +7t/2- or — 7r/2-flux on all triangles Ti.'^ In the 
following we will fix this flux to be = +7r/2. On the tri- 
angles and which involve the hopping, we 



FIG. 4. (a) A vortex (blue squares) superlattice of spacing D (here 
D = 2) and (b) the corresponding effective model in terms of Majo- 
rana zero modes tunneling on the triangular superlattice. The effec- 
tive hopping amplitudes ti (solid lines) and (dashed lines) are 
identified with the energy splittings e{D) and e(D\/3) at the corre- 
sponding vortex separations, respectively. The fluxes in the effec- 
tive model are chosen such that plaquettes indicated as Ti, T^, and 
Tj^ ^ are assigned fluxes <1>Ti = f , ^T_y^ — ~f ' ^T-^ ^/s ~ f ' 
respectively, that are consistent with the enclosed area. 

fix the fluxes to be = 37r/2 = —tt/2 (mod 2tt) and 
<i>Tj ^ = +7r/2 such that the flux through each triangle is 
proportional to its enclosed area. A periodic pattern of gauge 
choices Sij satisfying these flux assignments requires a mag- 
netic unit cell of at least 36 sites (see Appendix A). 

We (numerically) diagonalize this effective model for vary- 
ing relative amplitudes of the two hopping terms ti and t^. 
Fig. |5] illustrates the resulting phase diagram that shows var- 
ious gapped phases characterized by Chem numbers j^m = 
±1. ±3. The occurrence of some of these phases can be read- 
ily understood in the following way: If the nearest-neighbor 
hopping dominates {ti ^ t^)^ we recover the analytically 
tractable triangular lattice problem studied before'^, which 
has been shown to give rise to a gapped phase with = ±1. 
Similarly, in the opposite limit {t^ ^ ti) the system simply 
decomposes into three uncoupled copies of the triangular lat- 
tice problem, thus giving a Chem number of — ±3. There 
are also intermediate phases with Chern number = ±3 
that arise for —2 < ti/t^ < — 1 and arise due to a competi- 
tion between the two limiting cases above. The general sym- 
metry of the phase diagram where I'm —vm as t; -> — is 
characteristic for the triangular lattice, where the elementary 
plaquettes are odd cycles. Inverting the signs of all tunnel- 
ing amplitudes is equivalent to inverting flux ±7r/2 — >■ 
on all plaquettes, which leads to a time reversed phase with a 
Chern number of opposite sign. 

Equipped with these quantitative results for the effective 
Majorana model, we return to establishing a direct connec- 
tion between the microscopic pairwise vortex interactions and 
the collective ground state in the presence of a spacing D su- 
perlattice of interacting vortices. The interactions enter ([3]l by 
identifying the hopping amplitudes ti and with the corre- 
sponding energy splittings e{D) and e{D^/?>), respectively, as 
illustrated in Fig. [4] More precisely, we will use the ansatz 

ti = {-lf^-\e{lD)\, (4) 

where Pd is the fermionic parity {Pd = for even, Pd = 1 for 
odd) of the respective two vortex sector with vortex separa- 
tion d. The fermionic parity determining the sign of the phys- 



coupling parameters -h/i^^ 

FIG. 5. Left: The phase diagram of the effective Majorana 
fermion zero mode model (here ti = cos^ and — sin 6), as 
characterized by the Chern numbers for the flux assignment 
($Ti , $T^, $i,T^) = (f,— f,f). The squares and the circles 
denote phase transitions at jtij = \t^\ and ti — —2t^, respec- 
tively. Right: The lowest gap in the energy spectrum. Am, along the 
path shown in (a) for = — 1. 

ical energy splitting originates from the mapping of the spin 
Hamiltonian ([TJ to Majorana Hamiltonian on the honeycomb 
lattice'** and it is thus specific only to the honeycomb model. 
In the other microscopic models with interacting anyons both 
the magnitude and the sign can be obtained directly from the 
oscillating energy splitting. 

V. RESULTS 

To show that our effective model correctly predicts the 
Chern numbers of the nucleated phases, we will restrict to 
isotropic = Jy — couplings where the Hamilto- 
nian ([T]i has rotational C3 symmetry. Anisotropic couplings 
lead to anisotropic interactions, which require a more compli- 
cated effective model with correspondingly anisotropic hop- 
ping amplitudes.-" Thus, for simplicity we restrict ourselves 
to tuning only the magnitude K of the three spin term. It 
respects the rotational symmetry and the isotropic effective 
model (|3]l is valid for all values of it. The idea is then as 
follows: By tuning the microscopic coupling K we can tune 
the honeycomb model in the presence of a vortex superlattice 
through various topological phases. Independently this will 
tune also the pairwise energy splittings e{d) to different mi- 
croscopic values, which in turn through (j4]) will tune the state 
of the effective Majorana model. If the interactions are indeed 
responsible for the emergent Abelian phases, we expect agree- 
ment between the observed (i/) and predicted (C/m = vm — 1) 
Chern numbers. Furthermore, we expect the interaction in- 
duced nucleated gap (Am) to approximate the observed en- 
ergy gap (A^,). 

Our results are summarized in Fig. |6j which indeed shows 
quantitative agreement between the observed and predicted 
energy gaps and Chern numbers. This is most clearly il- 
lustrated for the D — ?> superlattice, where the observed 
sequence of phases with Chern numbers v = —4, 2, —2 
matches exactly that predicted by the effective model. In gen- 
eral, the approximation provided by the effective model gets 
more accurate for larger K and sparser vortex superlattices 
(for D — 1 the Chern numbers agree only for large K, while 
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FIG. 6. The energy gaps A„(A') of the nucleated 

phases in the honeycomb model (black circles) and the gaps 
/!^M\e{D,K),e{D\/i,K)] predicted by the effective model (red 
squares) for D = 1,2 and 3 vortex superlattices. is the Chem 
number observed in the honeycomb model and Dm = vm — 1 is 
the one predicted by the Majorana model from the interactions. Note 
the quantitative similarity of the sequence of phases for D = 3 su- 
perlattice to that of Fig. |5] (the finite gap around K — 0.0585 is a 
finite-size effect). All the data is for J^, = Jy = Jz = 1/3 and 
the values of the pairwise energy splittings e(d) have been obtained 
by restricting to the relevant two vortex sectors in a finite system of 
40 X 40 plaquettes (3200 sites) on a torus. 



for D > 2 they always agree). The reason is that for small K 
or for tight superlattices the coherence length ^ of the under- 
lying non-Abelian phase can increase (as the gap A ~' K^^ 
decreases) beyond the superlattice spacing^. This presumably 
renders the notion of individual vortices poorly defined and 
thus makes our microscopic approach inapplicable (while nu- 
cleation still occurs, the tunneling amplitudes are no longer 
captured by Q and/or many-body effects become relevant). 
Indeed, even if the nucleated gap A^, decays exponentially 
with increasing D, in agreement with it being induced by the 
interactions, further studies for spacings up to Z) = 6 shows 
excellent agreement over a wide range in K. 

Finally, we can now understand that the phases with Chern 
number v = —4 and 2 emerge due to the longer range inter- 
actions. Even if they are exponentially suppressed, the os- 
cillations can cause nearest neighbor interactions to vanish 
and thus make the next nearest interactions to dominate the 
physics. We could have included in the model interactions of 
even longer range, but as ([3]) can account for all the observed 
Chern numbers, we regard it as providing a complete descrip- 



tion. The exponential decay of the interactions means that if 
stable phases due to even longer range interactions existed, 
they should have appeared for the considered dense superlat- 
tices. 

VI. CONCLUSIONS 

We have studied in microscopic detail the nucleation of 
topological liquids - a general mechanism where interactions 
between non-Abelian anyons cause a transition to another 
topological phase^ - in the context of Kitaev's honeycomb 
model. Generally, we find that the presence of a superlattice 
of interacting Ising vortices always destroys the non-Abelian 
nature of the underlying topological phase and nucleates an 
Abelian topological phase. Employing an effective model of 
Majorana zero modes tunneling on the vortex lattive, we show 
that the character of the collective many-vortex state and the 
energy gap protecting it can be fully traced back to the sign, 
amplitude and ratio of the pairwise vortex-vortex interactions 
of different range. This provides an explicit demonstration 
of how local microscopies can change the global topological 
nature of the system. 

We found Abelian phases with Chem numbers ly — +2 
or V = —4 that arise when oscillations in the interactions, 
even if they decay exponentially, cause next-nearest neigh- 
bor interactions to become of comparable (but of different 
sign) or larger magnitude than the nearest-neighbor interac- 
tions. Due to the omnipresence of such oscillating anyon- 
anyon interactions,^^ our results imply that a similar nu- 
cleation mechanism should occur also in other microscop- 
ically distinct topological systems. These include Moore- 
Read quantum Hall states, where quasihole excitations can 
form a Wigner crystal when the magnetic field is detuned 
away from the middle of the plateau,'' p-wave superconduc- 
tors, where vortices can form an Abrikosov lattice, and topo- 
logical nanowires arranged in regular arrays.'^ While we fo- 
cused here on nucleation in the presence of an ideal vortex 
superlattice, further studies^" will address the microscopies 
of a recently observed thermal metal state in the presence of 
disorder.— 
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Appendix A: Unit cell for the Majorana model with longer 
range tunneling 

When = the unit cell with <^Ti = 7r/2-flux on all 
plaquettes contains two sites. In the other limiting case when 
ti = 0, the unit cell with = — 7r/2 consists of six sites, 
two each from each of three disjoint sublattices on which the 
next nearest hopping act. When both are present, one needs 
also to fix the flux on all the intermediate plaquettes, that con- 
sist of two ii -links and a single i^g-link, to be ^t-^ ^ = 7r/2. 
It turns out that a gauge consistent with all these flux assign- 
ments is translationally invariant only with respect to a mag- 
netic unit cell of 36 sites, as illustrated in Fig|7] 

For each site in the unit cell one has six independent 
gauge choices - {sf)ij,{s\)ij and (sf )ij on the ti-lattice and 



(s^)ij,('Sy3)ij and (s^)ij on the t^g-lattice. Taking into 
account the orientation of the links when evaluating the flux 
(see Fig|7]), one possible gauge choice giving rise to the de- 
sired (<I>Ti,$T^,$i,T^) = (7r/2, -7r/2, 7r/2) flux pattern 
is given by 
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FIG. 7. Left: The unit cell for the effective Majorana model with 
both nearest (ii) and next nearest (iyj) neighbor tunneling. The first 
couple all sites of the lattice, while the latter couple only sites belong- 
ing to one of three distinct sublattices denoted by circles, squares and 
diamonds. When — the unit cell consists of a black and a white 
site iiTespective of their sublattice label, while for ti = the unit cell 
consists of six sites - a black and a white site from each sublattice. 
Right: For every site {i,j) in the unit cell one has six independent 
gauge choices - three on the ti lattice and three on the respective 
sublattice. When calculating the flux per plaquette, one has to 
take into account the overall i factor in the tunneling. We assume 
a convention that the phase of the hopping between any two sites is 
given by i (— i) when it is along (against) the shown orientations. 



